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We construct a variety of novel localized states with distinct topological structures in the 3D 
discrete nonlinear Schrodinger equation. The states can be created in Bose-Einstein condensates 
trapped in strong optical lattices, and crystals built of microresonators. These new structures, 
most of which have no counterparts in lower dimensions, range from purely real patterns of dipole, 
quadrupole and octupole types to vortex solutions, such as "diagonal" and "oblique" vortices, with 
axes oriented along the respective directions (1,1,1) and (1,1,0). Vortex "cubes" (stacks of two 
quasi-planar vortices with like or opposite polarities) and "diamonds" (discrete skyrmions formed 
by two vortices with orthogonal axes) are constructed too. We identify stability regions of these 
3D solutions and compare them with their 2D counterparts, if any. An explanation for the stabil- 
ity/instability of most solutions is proposed. The evolution of unstable states is studied as well. 



Introduction. In the past two decades an explosion of 
activity has been observed in the study of intrinsic lo- 
calized modes (ILMs, alias discrete solitons) in nonlinear 
dynamical lattices, especially due to the ability of such 
modes to act as energy "hot spots" lj. The relevance 
of ILMs has been demonstrated in problems ranging 
from arrays of nonlinear-optical waveguides and pho- 
tonic crystals |3j, to Bose-Einstein condensates (BECs) 
trapped in optical lattices (OLs) [j] and Josephson- 
j unction ladders [5(. 

A universal model, which may arise as an envelope ap- 
proximation from most of the complex nonlinear equa- 
tions on the lattice and also as direct physical model for 
the BECs 4] and optical waveguiding arrays [a, IS, is the 
discrete nonlinear Schrodinger (DNLS) equation [g. On 
top of its significance to applications, the DNLS equation 
itself is a fundamentally interesting dynamical model. In 
the 3D case, its direct physical realization is provided, 
as mentioned above, by BECs trapped in strong OLs J4j- 
Waveguide arrays, however, cannot be described by a 3D 
discrete model, since the evolution variable in the op- 
tical media is a spatial coordinate, while the temporal 
variable, which effectively plays the role of an additional 
quasi-spatial one, cannot be discrete. Nevertheless, an- 
other physical realization of the 3D DNLS equation may 
be provided by a crystal built of microresonators . 

The study of the 3D continuum NLS equation, includ- 
ing a 3D Jm or quasi-2D 11] OL, and of the DNLS model 
proper |12|. has started recently, becoming increasingly 
accessible to numerical computations. As a result, the 
first coherent structures, such as discrete vortices of the 
topological charge (vorticity) S =1,2 and 3, were identi- 
fied and their stability was investigated. The aim of the 
present work is to study a large variety of novel localized 
3D structures in the DNLS equation, many of which turn 
out to be stable. In particular, we construct states which 
include, first, dipoles with the axis oriented along a lattice 
bond, or along a planar diagonal, or along a 3D diago- 



nal (we call them, respectively, "straight" , "oblique" , and 
"diagonal" dipoles). Next, we construct quadrupole and 
octupole states, that, similarly to the dipoles, are real so- 
lutions. More sophisticated structures are also presented, 
namely "vortex cubes" (concatenations of two straight 
vortices with the same or opposite charges centered on 
parallel planes), oblique and diagonal vortices, and "vor- 
tex diamonds" , formed by a crossed pair of vortices with 
orthogonal axes. Apart from the straight and oblique 
dipoles and quadrupoles, these ILMs have no counter- 
parts in 2D lattices. 

To present the results, we first introduce the model, 
and then report systematic numerical results for the 
shape and stability of the new localized states. This is 
followed by conclusions, which include an explanation for 
the stability and instability of the majority of patterns 
found in this work. 

The Model. We consider the DNLS equation on the 
cubic lattice with a coupling constant C |12j ] , 

o, (i) 

with (j) = d(f>/dt, and the discrete Laplacian is A(j)^ m ^ n = 

<Pl+l,m,n + 4>l,m+\,n + <f>Lm,n+l + 4>l-l,m,n + 4>l,m-l,n + 

4>i,m,n-i — §4>i,m,n ■ Solutions are looked for as 4>i t m,n = 
exp(iAi) with a frequency -A (or the chemical po- 
tential in the context of BEC), where the stationary func- 
tions ui jm ^ n obey the equation 

Alt; !?Tl , n = CAu/ im! „ + |li; !m . n | Ui^ m ^ n . (2) 

Profiles used as an initial guess for the fixed-point itera- 
tion converging to solutions displayed below, were based 
on the form of the respective solutions (for the same A) 
in the anti-continuum (AC) limit, C = 0. Once solutions 
to Eq. have been obtained, the linear-stability anal- 
ysis is performed for a perturbed solution [12j, (f>i, m ,n = 
[ui,m,n + e (ai, m , n e Xt + ii, m , n e A *')] e lAt , where e is an in- 
finitesimal amplitude of the perturbation, and A is its 



eigenvalue. The Hamiltonian nature of the system dic- 
tates that if A is an eigenvalue, then so are —A, A* and 
—A* (in the stable case, A is imaginary, hence this sym- 
metry yields only two different eigenvalues, A and —A). 
The stationary solution is unstable if at least one pair of 
the eigenvalues features nonvanishing real parts. 

Results. We start by constructing purely real solutions 
of the dipole type (in lower-dimensional models, solitons 
of this type were considered in Refs. 0, Figure 
n displays generic examples of "tight" dipoles with ad- 
jacent excited sites (i.e., the separation between them 
is d = 1, in the corresponding units) and three possi- 
ble orientations relative to the lattice: straight (along a 
lattice's bond), (a), oblique (along a planar diagonal), 
(b), and diagonal (along a 3D diagonal), (c). In this fig- 
ure and below, unless stated otherwise, we show a typ- 
ical case admitting stable solutions, with C — 0.1 and 
A = 2. The borders of the stability windows, < C < 
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, for the three types (straight, oblique and diago- 



nal) of the dipoles are given by C^^ tl . = 0.23013 ± SC, 

Cff-cL = 0.53666 ± <5C7, and = 0.73084 ± SC, 

where the error margin is SC = 0.00001 for all the results 
presented in this work (unless indicated otherwise). We 
observe that the diagonal dipole in Fig. [TJc) remains sta- 
ble in a larger interval than its oblique counterparts in 
Fig. [Hb) that, in turn, is more stable that the diagonal 
one in Fig. ^c). It is also possible to construct dipole 
solutions corresponding to d = 2, with excited sites sepa- 
rated by a single nearly empty site. Such dipole solutions 
(not depicted here) have larger stability windows than 



their "tight" (d = 1) counterparts: C di 



(3D, 2) 



0.66722, 



°dip-obl 



p— str 

1.08024, and = 1.31356. As seen 

in Fig. nih), C^fP'^ further increases with the sepa- 
ration distance a, approaching the stability threshold 
of the fundamental (single-site-based) discrete soliton, 
c £sd,oo) ^ c jSD) = 2 oog ± Q Q01 jj^ gee the horizontal 

dashed line in the figure. It can also be found that the 

same relations < T (3D ' d) < T (3D ' d) < r (3D4) as 
same relations, U < o dip _ str < o dip _ obl < o dip _ dia , as 

found for d — 1 , are valid for all d, see Fig. njh) . An- 
other relevant comparison is with dipoles in the 2D DNLS 
model, which were studied in Rcf. hj. The comparison 
shows that the dipole soliton of the oblique and straight 
types, which have their counterparts in the 2D case, are 
less stable, although not drastically, than those counter- 
parts: cgj^ = 0.23013 < tfg^ = 0.245 ± 0.005 
for A = 2. The weaker stability of the 3D structures is 
explained by the analogy with the multidimensional con- 
tinuum NLS equation, where all the solitons are destabi- 
lized by collapse, which is, respectively, weak and strong 
in the 2D and 3D cases [l6j. 

Quadrupole and octupole solitons are also shown in 
Fig. The quadrupole is based on four contiguous sites 
(so that we prescribe d = 1 to this structure too) which 
form a square in the plane, Fig. IHd). It was found to 

be stable for C < C^ ] = 0.13836, while its 2D analog 
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FIG. 1: (Color online) Stable multipoles. The top row de- 
picts stable "tight" dipoles (with d — 1): (a) straight, (b) 
oblique, and (c) diagonal ones, (d) and (e): Quadmpoles in 
the n = plane, with the internal separation d = 1 and d = 2, 
respectively, (f) and (g): Octupoles with d = 1 and d = 2. 
Panel (h) displays the stability threshold cir U ' as a function 
of the internal distance d for (from top to bottom) diagonal, 
oblique, and straight dipoles, octupoles and quadrupoles. The 
horizontal dashed line corresponds to the stability threshold 



for the fundamental discrete soliton, C, 



(3D) 



2. Note that, 



for the quadrupole (the bottom graph), C^jf^ behaves lin- 
early for small d (see the dashed line with slope 0.325 for 
guidance). In panels (a)-(g), level contour corresponding to 
Re(ui,m,n) = ±0.5 are shown in blue and red (dark gray and 
gray, in the black-and-white version), respectively. All these 
states are stable (for the case shown, with A = 2 and C — 0.1). 



has C*^ d 1} = 0.1485 ±0.0005 Q . In this case, as well as 
it was with the dipoles, the 2D configurations tend to be 
slightly more stable than their 3D siblings. The octupole 
is shown in panel ^f); it is based on a set of eight con- 
tinuous sites (therefore it is also assigned d = 1) forming 
a cubic cell in the 3D lattice. It is stable in the interval 
C < C^± = 0.10030, which is smaller than the above 
ones for the quadrupoles and dipoles. Similar to the 
dipoles, multipoles can also "swell" by inserting unpop- 
ulated sites between the excited ones. The resulting sta- 
bility intervals for d — 2 are larger than their d — 1 coun- 
terparts: C^f = 0.503232 and C 3 °' 2) = 0.418411, see 
Figs.QIe) and (g), respectively. It is relevant to mention 
that all the newly found structures have their stability 
limits lower than those for the fundamental (single-site- 
based) discrete soliton, C^ = 2.009 ± 0.001 [13, see 
Fig. Eh). 

Stable dipoles have one pair of imaginary eigenvalues 
in their perturbation-mode spectrum, that, with the in- 
crease of C, collides with the continuous spectrum, which 
leads to the destabilization. On the other hand, stable 
quadrupoles have three such pairs (two of them form a 



FIG. 2: (Color online) Vortex cubes for A = 2 and C = 0.1. 
Panel (a) shows a stable vortex cube, built of two quasi-planar 
vortices, set in the planes n = ±1, with equal vorticities, Si = 
£2 = 1, and a phase shift of ir. Panel (b) shows an unstable 
cube formed by vortices with opposite charges, Si = —52 
and panel (c) shows a snapshot (at t = 200) of its evolution, 
clearly demonstrating that the phase coherence between sites 
forming the pattern is lost. The real level contours are as 
in Fig. and the imaginary ones, Im(ui t m,n) = ±0.5, are 
shown by green and yellow (light and very light gray, in the 
black-and-white version) hues, respectively. 

doublet for small C) , and the octupoles have seven (six of 
which form two triplets for small C) . More generally, the 
number of potentially unstable eigenvalue pairs is N — 1, 
where N is the number of sites on which the structure is 
based 0. 

The next novel type of a 3D discrete soliton, with no 
lower-dimensional counterpart, is a "vortex cube", which 
is built as a stack of two quasi-planar vortices with equal 
topological charges S\ = S2 = 1 and a phase shift A(j> = 
it, separated by an empty layer, so that it has d = 2. 
Figure|2Ia) shows real and imaginary parts of the vortex- 
cube lattice field. Such a state is stable for < C < 
^cub 7^ = 0-56324. On the contrary to this, a vortex cube 
built as a stack of two in-phase vortices, with A<fi = 0, 
is always unstable, through three real eigenvalue pairs. 
Further, Fig.^b) shows a similar stack, but composed of 
two vortices with opposite charges, S\ = — S2 = 1. This 
configuration is always unstable as well, due to a real 
eigenvalue pair. In this case, the instability manifests 
itself as a symmetry breaking between the two planes 
and, hence, the phase coherence of the entire pattern is 
eventually lost, see Fig.^c). 

Another 3D object, with no lower-dimensional analog 
either, is a vortex with the axis directed along the di- 
agonal of the cubic lattice, i.e., the vector (1, 1, 1). Fig- 
ure|2a) shows such a "diagonal vortex" constructed by a 
continuation procedure starting, in the AC limit, with the 
following distribution of local phases: 4>i,m,n = litSk/Q 
(k = 0, 1, 2, 5), for the sites that lie in a plane orthogo- 
nal to the axial (diagonal) direction: (1,-1, 0), (0, —1, 1), 
(-1,0,1), (-1,1,0), (0,1,-1), (1,0,-1). Obviously, this 
phase pattern bears the vorticity S (S = 1 in Fig. |3J). 
However, the diagonal vortex turns out to be always un- 
stable, due to three real eigenvalue pairs, and it eventually 
settles to a single-site-based soliton. 

One more species of discrete vortices that may ex- 
ist solely in the 3D lattice is an "oblique" one, shown 
in Figs. EJb)-(c), with the axis' directed along (1,1,0). 
In the AC limit, the solution is carried by the array of 
sites (1,-1,0), (1,-1,1), (0,0,1), (-1,1,1), (-1,1,0), 



FIG. 3: (Color online) Diagonal and oblique vortices. Panel 
(a) shows an unstable "diagonal" vortex, with the axis along 
the direction (1, 1, 1), for A = 2 and C = 0.1. Panel (b) shows 
an unstable oblique vortex with the axis oriented along the 
direction (1, 1,0), and panel (c) shows a stable oblique vortex 
of a modified form (see text) for A = 2 and C = 0.01. 




FIG. 4: (Color online) The vortex diamond (discrete 
skyrmion) for A = 2 and C = 0.1. Panel (a) displays an 
unstable "diamond", and panel (b) depicts its evolution, in 
terms of the field's magnitude (top) and phase (bottom) at 
the main six lattice sites. Initially [t < 300), pairs of sites (op- 
posite vertices of the diamond) cyclically change the phase, 
and subsequently (t > 350) the phase coherence is lost and 
the field magnitude oscillates erratically about its initial level. 

(—1,1,-1), (0,0,-1), and (1,-1,-1), with the phases 
at them S ■ (0, a, tt/2, a + tt/2, tt, a + tt, 37r/2, a + 3tt/2), 
where a = tan _1 (l/v / 2). Figure |3fb) depicts such an 
oblique vortex, which, in this form, is found to be always 
unstable, similar to the diagonal vortex. Nonetheless, the 
oblique vortex can be stabilized in a modified form, by in- 
troducing a sign shift at the intermediate edge sites [i.e., 
(0, 0,1), (-1,1,0), (0, 0,-1) and (1,-1,0)], see Fig.^c). 
This sign change avoids contiguous sites with the same 
phase and does not alter the vortex' topological charge, 
which remains 1. The modified oblique vortex is stable 
in a small interval, C < C^f2 ob] = 0.0104. 

Finally, motivated by the concept of skyrmions |19| . 
which are distinguished by two topological charges asso- 
ciated with closed contours in two perpendicular planes, 
we have constructed one more type of vortex structures 
in the 3D lattice, viz., a "diamond" shown in Fig.^Ja). 
It is built as a vortex cross, i.e., a nonlinear superposi- 
tion of two straight site-centered S = 1 vortices, with 
their axes directed along two orthogonal directions. The 
stability analysis shows that, although the diamond has 
three imaginary (stable) eigenvalue pairs, it is always un- 
stable due to a real eigenvalue pair. In direct simula- 
tions, its instability manifests itself in a rather intriguing 
manner, see Fig. 0Jb). At t 1=3 100, two pairs of oppo- 
site vertices of the diamond change phases (one by +ir/2 
and the other by —tt/2) and suffer a momentary mag- 
nitude change. After this, the diamond remains stable 
until t « 200, when the same pairs of sites suffer an- 
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other phase shift (in the same direction). This process 
repeats itself almost periodically until (at t « 340) the 
diamond, after 4 shifts, returns back to its original phase 
distribution. Subsequently, the phase shifts accelerate 
and phase coherence is finally lost. Simultaneously, am- 
plitude variations of ±10% get accumulated, see the top 
panel of Fig. Efb). Eventually, the solution degenerates 
into a plain single-site-based soliton. 

Conclusions and discussion. We have introduced sev- 
eral novel species of topologically structured discrete soli- 
tons in 3D dynamical lattices, using the paradigm of the 
discrete nonlinear Schrodinger equation. The solutions 
have been constructed starting from properly chosen 
anti-continuum approximations, and their linear stabil- 
ity was studied through the computation of the relevant 
eigenvalues. Previously, only the fundamental single-site- 
based solitons and straight vortices, with the axis di- 
rected along a lattice bond, were known. We have found 
three species of dipoles, which differ by the orientation 
relative to the lattice, quadrupoles and octupoles, vor- 
tex cubes (stacked dual- vortex patterns), diagonal and 
oblique vortices, and "diamonds" (vortex crosses or dis- 
crete skyrmions). Except for the straight and oblique 
dipoles and quadrupoles, the patterns obtained are en- 
demic to the 3D lattice setting, having no counterparts 
in lower dimensions. 

Apart from the diagonal vortices and "diamonds" , all 
the patterns constructed above have stability regions be- 
low a critical value of the coupling parameter. It is pos- 



sible to explain the stability/instability of all the struc- 
tures that may be realized as bound states of two simpler 
objects, viz., dipoles, quadrupoles (bound states of two 
dipoles with opposite orientations), octupoles (bound 
states of two quadrupoles), and vortex cubes. Indeed, a 
known general principle is that a bound state pinned by 
the lattice may be stable only if the coupled objects repel 
each other [l5|,|20( (i.e., have a phase difference of ir be- 
tween their building blocks). This explains the existence 
of stability regions for multipoles of all types. Similarly, 
considering the interaction between constituent quasi- 
planar vortices, one may understand the stability and 
instability of vortex cubes of the types shown in Figs. 
G^a) and E{b), respectively. Following this principle, it is 
also possible to predict the stability of more exotic 3D 
patterns such as bound states of two oblique or diago- 
nal dipoles, or octupoles constructed of two such states. 
In those cases when the 3D structures have 2D counter- 
parts, viz., straight and oblique dipoles and quadrupoles, 
their stability regions are narrower than in the 2D case, 
which is explained by a stronger trend to collapse in three 
dimensions. 

Future challenges involve semi-analytical investigation 
of such solutions via Lyapunov-Schmidt theory, and iden- 
tification of their stability by means of methods similar 
to those developed for the ID and 2D cases [l^. On the 
experimental side, it may be interesting to create such 
structures in BEC loaded into a strong optical lattice. 



[1] S. Aubry, Physica 103D, 201 (1997); S. Flach and C.R. 

Willis, Phys. Rep. 295, 181 (1998); D.N. Christodoulides 

F. Lederer, and Y. Silberberg, Nature 424, 817 (2003); 

D.K. Campbell, S. Flach, and Y.S. Kivshar, Phys. Today, 

January (2004), p. 43. 
[2] see e.g., A. A. Sukhorukov, Y.S. Kivshar, H.S. Eisenberg, 

and Y. Silberberg, IEEE J. Quant. Elect. 39, 31 (2003); 

U. Peschel et al, J. Opt. Soc. Am. B 19, 2637 (2002). 
[3] S. F. Mingaleev and Y.S. Kivshar, Phys. Rev. Lett. 86, 

5474 (2001). 

[4] A. Trombettoni and A. Smerzi, Phys. Rev. Lett. 86, 2353 
(2001); F.Kh. Abdullaev et al., Phys. Rev. A64, 043606 
(2001); F.S. Cataliotti et al, Science 293 , 843 (2001); 
A. Smerzi et al, Phys. Rev. Lett. 89, 170402 (2002); G.L. 
Alfimov, P.G. Kevrekidis, V.V. Konotop, and M. Salerno, 
Phys. Rev. E 66, 046608 (2002); R. Carretero-Gonzalez 
and K. Promislow, Phys. Rev. A 66 033610 (2002). 

[5] P. Binder et al, Phys. Rev. Lett. 84, 745 (2000); E. Trias 
et al, Phys. Rev. Lett. 84, 741 (2000). 

[6] D.N. Christodoulides and R.I. Joseph, Opt. Lett. 13, 794 
(1988); A. Aceves et al, Opt. Lett. 19, 1186 (1994). 

[7] H.S. Eisenberg et al, J. Opt. Soc. Am. B 19, 2938 (2002); 
J. Meier et al, Phys. Rev. Lett. 91, 143907 (2003). 

[8] P.G. Kevrekidis, K.0. Rasmussen and A.R. Bishop, Int. 
J. Mod. Phys. B 15, 2833 (2001). 



[9] J.E. Heebner and R.W. Boyd, J. Mod. Opt. 49, 2629 

(2002); P. Chak et al, Opt. Lett. 28, 1966 (2003). 
[10] B.B. Baizakov, V.V. Konotop, and M. Salerno, J. Phys. 

B 35, 5105 (2002); B.B. Baizakov, B.A. Malomed, and 

M. Salerno, Europhys. Lett. 63, 642 (2003). 
[11] B.B. Baizakov, B.A. Malomed, and M. Salerno, Phys. 

Rev. A 70, 053613 (2004). 
[12] P.G. Kevrekidis, B. Malomed, D.J. Frantzeskakis and R. 

Carretero-Gonzalez, Phys. Rev. Lett. 93, 080403 (2004). 
[13] P.G. Kevrekidis, B.A. Malomed, A.R. Bishop, and D.J. 

Frantzeskakis, Phys. Rev. E 65, 016605 (2002). 
[14] J. Yang et al, Opt. Lett. 29, 1662 (2004). 
[15] P.G. Kevrekidis, B.A. Malomed, and A.R. Bishop, J. 

Phys. A 34, 9615 (2001). 
[16] C. Sulem and P.L. Sulem, The Nonlinear Schrodinger 

Equation, Springer- Verlag (New York, 1999). 
[17] P.G. Kevrekidis, B.A. Malomed, Z. Chen, and D.J. 

Frantzeskakis, Phys. Rev. E 70, 056612 (2004). 
[18] D.E. Pelinovsky, P.G. Kevrekidis, and D.J. Frantzeskakis 

nlin.PS/0411016 and nlin.PS/0410005 
[19] T.H.R. Skyrme, Proc. Roy. Soc. London A 260, 127 

(1961). 

[20] T. Kapitula, P.G. Kevrekidis, and B.A. Malomed, Phys. 
Rev. E 63, 036604 (2001). 



